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Abstract 

We presents results for the matrix elements of a variety of four-fermion operators 
calculated using quenched Wilson fermions. Our simulations are done on 170 lattices 
of size 32 3 x 64 at (3 = 6.0. We find B K = 0.74 ± 0.04 ± 0.05, B D = 0.78 ± 0.01, 
£? /2 = 0.58±0.02+J;JJ, Bl' 2 = 0.81±0.03 + ° 3 2 , with all results being in the NDR scheme 
at ji = 2 GeV. We also calculate the B-parameter for the operator Q s , which is needed 
in the study of the difference of B-meson lifetimes. Our best estimate is Bg(NDR, /i = 
1/a — 2.33 GeV) = 0.80 ± 0.01. This is given at the lattice scale since the required 2- 
loop anomalous dimension matrix is not known. In all these estimates, the first error is 
statistical, while the second is due to the use of truncated perturbation theory to match 
continuum and lattice operators. Errors due to quenching and lattice discretization are 
not included. We also present new results for the perturbative matching coefficients, 
extending the calculation to all Lorentz scalar four-fermion operators, and using NDR 
as the continuum scheme. 
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1 Introduction 



One of the central goals of lattice QCD is to calculate hadronic matrix elements of phenomenological 
interest [l], § . We present here results for the matrix elements of a variety of four-fermion operators: 
Bk-i which is needed as input for estimates of CP violation in K—K mixing; Bp, required to estimate 
D — D mixing, and an indication of the result for B mesons; Bj~ 3 ^ 2 and B^~ 3 ^ 2 , which determine 
the dominant contribution of electromagnetic penguins to e'/e; and a new quantity, B$, needed as 
part of an estimate of the lifetime difference between B mesons 0. 

Our results are obtained using Wilson fermions in the quenched approximation. The use of Wil- 
son fermions makes it more difficult to extract quantities which are constrained by chiral symmetry. 
This is a problem for the calculation of B^, and our results are not competitive with those from 
staggered fermions. Our aim here is to gain experience at reducing systematic errors. The other 
matrix elements are not constrained by chiral symmetry, and for these it is more straightforward 
to use Wilson than staggered fermions. 

Compared to previous work with Wilson fermions, our study uses larger lattices (32 3 x 64), 
which should make finite volume effects negligible and reduces the statistical errors, has improved 
statistics (170 lattices), and uses a larger range of light quark masses (0.3m S)P hys-2m SiP hys)- The 
latter improvement allows us to do more reliable chiral extrapolations, since we can include terms 
of 0{m A K ). We have results at only a single lattice spacing, (3 = 6/g 2 = 6, so we cannot extrapolate 
to the continuum limit of the quenched theory. What we provide is one point with statistical errors 
small enough that systematic errors due to chiral extrapolations and the truncation of perturbation 
theory can be quantified. 

As an adjunct to our numerical results, we have extended previous calculations of the perturba- 
tive matching coefficients between lattice and continuum operators. We also present a renormaliza- 
tion group improved matching formula which incorporates tadpole improved lattice perturbation 
theory. 

The plan of the paper is as follows. We begin, in Sec. 0, by summarizing our results. In Sec. |3] 
we describe the various technical details pertinent to this calculation. Section |] summarizes the 
matching between lattice and continuum operators. Results for the various ^-parameters are given 
in the final four sections. We discuss two technical issues in appendices. Appendix A describes a 
renormalization group improved matching formula, while in Appendix B we present general formulae 
for the one-loop matching coefficients for all Lorentz scalar four-fermion operators between the 
lattice and the NDR scheme. 

2 Summary of results 

We begin with our result for B K (defined in Eq. 

B K (NDR, 2 GeV) = 0.74 ± 0.04(stat) ± 0.05(subt) . (2.1) 

The second error is an estimate of the uncertainty due to our method of subtracting chiral artifacts. 
This is a considerable improvement over our previous work [|J — the statistical error has been reduced 
by a factor of ~ 6, and we have better control over the systematic errors. 

It is interesting to compare our result to those other recent high statistics calculations at or 
near the same lattice spacing, particularly since all use different subtraction methods. Bernard and 
Soni use Wilson fermions at (3 = 6 on 24 3 x 39 lattices 0, and find B^ = 0.67 ±0.07. Their 
subtraction method uses off-shell matrix elements |J whereas ours uses on-shell matrix elements 
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at finite momentum transfer. One concern with their result is that they use larger quark masses, 
which, according to our analysis, may lead to an overestimate of B^. Nevertheless, it is encouraging 
that the results agree. 

The APE group has used improved "clover" fermions and a non-perturbative determination 
of the matching coefficients 0, |§. Their non-perturbatively matched operator has the correct 
chiral behavior. They present their result (from j3 = 6 on 18 3 x 64 lattices) in the "regularization 
independent" scheme, Bk(RI, 2.02 GeV) = 0.62(11). The numerical value is unchanged upon 
conversion to NDR, and is consistent with our result. The two results can differ by terms of 0(a), 
but the errors are too large to resolve any difference. 

The JLQCD collaboration has recently presented results using Wilson fermions and a non- 
perturbative determination of matching coefficients 0. They find, at (3 = 5.9 and 6.1, that 
B K (NDR, 2 GeV) = 0.48 ± 0.05 and 0.70 ± 0.07, respectively. They also have a result for the 
value in the continuum limit, Bk{NDR,2 GeV, a = 0) = 0.59 ± 0.08. The latter two numbers 
are consistent with our result, while the former is not. Our methods seem to work comparably 
well — the final errors are similar based on similar samples and lattice sizes. It will be interesting 
to use their non-perturbative matching coefficients with our data to evaluate the accuracy of the 
procedure we use to remove the lattice artifacts in the chiral expansion. 

Finally, we compare our result with that from staggered fermions. We expect that the two should 
agree up to corrections of 0(a). At (3 = 6.0, different choices of discretized staggered operators 



give results in the range 0.68 — 0.71 [JTOL O, |12, 13|, with statistical errors of 0.01 or smaller. Our 



results are consistent with these, but our errors are too large to allow us to see the expected 0(a) 
differences. 

Our final results for 5-parameters for D mesons (defined in Eq. |b.l| ) are 

5 D (NDR,/i = 2 GeV) = 0.785 ±0.015, (2.2) 
5fl s (NDR,/i = 2 GeV) = 0.83 ± 0.01 , (2.3) 
B Ds /B D = 1.047 ±0.014. (2.4) 

The errors are much smaller than for Bk as we need make no corrections for mixing with wrong 
chirality operators. Thus, for this calculation, we are at a point where the remaining systematic 
errors (finite lattice spacing and quenching) are likely to be larger than the statistical errors. These 
results are consistent with those previously obtained (see, for example, the compilation in Refs. []5| 
and |15[), but have considerably smaller errors. 



Our final results for the S-parameters of the 1 = 3/2 parts of the operators <2 7 and Q 8 (defined 
in Eqs. [TT]-^ are 

£?7^ 2 (NDR, 2 GeV) = 0.58 ± 0.02(stat) + °'%pert) , (2.5) 

U .Uo 

Eg /2 (NDR,2 GeV) = 0.81 ± 0.03(stat) + °'° 3 (pert) . (2.6) 

.02 

The "perturbative error" reflects the dependence of the results on the choice of a s used in the 
matching of continuum and lattice operators, and is comparable to or larger than the statistical 
errors. The perturbative error could be removed by the use of non-perturbative matching coeffi- 
cients, but these are not yet available for Qj and Q$- Our results are consistent with those we found 
previously H, the apparent difference being due to the use of a different continuum regularization 
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scheme, a different final scale, and a different choice of a s in the matching of continuum and lattice 
operators. 

Our values are smaller than the numbers used by Lusignoli et a/.[[L| and Ciuchini et aZ. [IB| in 
their analyses of predictions for e'/e: they take B^ 2 = B^ 2 = 1.0 ±0.1. This difference is important 
because a smaller B s means a larger e'/e. 

The final quantities we consider are Bs = and the related parameter B§ [defined in Eqs. 



and 8.2]. What is of interest for phenomenology is the value of Bg for bs mesons. The closest we 



can come is the result for b quark with roughly the mass of the charm, which is 

S^(NDR,l/o) = 0.80±0.01(stat) , (2.7) 
5+(NDR, 1/a) = 0.94 ± 0.01 (stat) . (2.8) 

We have to quote these results at « = 1/a = 2.33 GeV because the two-loop anomalous dimension 
matrix needed to run to 2 GeV has not been calculated. For the same reason, we cannot estimate 
a "perturbative" error. 

There are no previous results for these 5-parameters using propagating heavy quarks. Recently, 
Gimemez and Martinelli have calculated Bs in the static limit, nib — > oo [IT]]. Their lattice result, 
also at (3 = 6.0 but without including any perturbative corrections, is Bg M = — (5/8)Bs = — (0.60± 
0.01 ± 0.03). While this is m 20% larger than our result, the above mentioned differences preclude 
any meaningful comparison. 



3 Technical Details 

We use a sample of 170 gauge configurations of size 32 3 x 64 generated at (3 = 6 in the quenched 
approximation. In the following we outline the method we have used to calculate matrix elements 
of four-fermion operators. Further details of our update and inversion algorithms, and of our 
determination of quark masses, can be found in Ref. [O . 



3.1 Quark propagators 

In our approach we need two kinds of quark propagators: one which allows for the creation of mesons 
with an explicit zero-momentum projection, and the other that allows good overlap with a range 
of lattice momenta. For the former we use wall sources (on a time slice fixed to Coulomb gauge); 
for the latter, gauge invariant Wuppertal sources. Both are calculated using periodic boundary 
conditions in all directions. We calculate propagators at five quark masses: k = 0.135 (C), 0.153 
(S), 0.155 (Ui), 0.1558 (U 2 ), and 0.1563 (U 3 ). These quark masses correspond to pseudoscalar 
mesons of mass 2835, 983, 690, 545 and 431 MeV respectively. Here, as in the following, we have 
used 1/a = 2.330(41) GeV, the scale we determined in Ref. |L4[] using M p . 

With these five flavors we construct 15 distinct "kaons", whose matrix elements we then study. 
The three C/j quarks allow us to extrapolate to the physical isospin symmetric light quark mass 
m = (m u + m,i)/2, whose value is determined by matching M n /M p to experiment. The physical 
value of strange quark mass, determined using M^/Mp, lies between S and Ui, and we use these 
two points to interpolate to it. It turns out that other ways of determining m s , e.g. matching to 
M^/Mr, lead to values differing by ~ 20%. To estimate the uncertainty that this introduces, we 
calculate kaon matrix elements in two ways: by interpolating to our standard m s , and by directly 
using the results for U2U3, which turns out to have almost exactly the physical kaon mass, albeit for 
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almost degenerate quarks. For Bk^B^ 2 , and B^ 2 , the difference between the two results turn out 
to be much smaller than other errors, so we do not quote a systematic error due to the uncertainty 
in setting m s . In case of B Ds and B$, the difference between using m s (M K ) and m s (M ( f,) is ~ 1% 
with m s {M^) giving the larger value. For reasons explained in [14], we believe that m s {M ( j ) ) gives 
a better estimate of m s , and we use this value consistently for our final results. To illustrate the 
details of our analysis we often use results from U2U3 as they are very close to the physical kaon 
mass. 

Finally, we take the physical charm mass to be C. With this choice the experimental values of 
Md, Md* and Md s lie between the static mass Mi (measured from the rate of exponential fall-off 
of the 2-point function) and the kinetic mass defined as M 2 = (d 2 E/dp 2 \ p=0 )~ 1 (b 



3.2 Extracting the four-fermion matrix elements 

We use the same method for calculating matrix elements of four-fermion operators as in our earlier 
work 0]. The method is illustrated in Fig. [l|. The initial pseudoscalar state is created by wall 
sources at time T = 0, and thus has p = (0,0,0). It propagates both forward and backwards in 
time. The four-fermion operator is inserted at a time in the range 1 — 31 (63 — 33 for the backward 
moving particle). The operator insertion is done at the 5 lowest lattice momenta, p = (0,0,0), 
(1, 0, 0), (1, 1, 0), (1, 1, 1) and (2, 0, 0), averaging over all possible permutations of the components 
of p. The pseudoscalar emerging from the four-fermion operator is then destroyed by an operator 
constructed using Wuppertal sources at T = 32. It is essential that this operator have a large 
coupling to kaons with all the above momenta, and we find that Wuppertal smearing does this well. 
In the end, we have two measurements, corresponding to the forward and backward propagation, 
and, furthermore, each of these is an average over a certain number of time slices. 

As shown in Fig |I], we calculate the ratio of the matrix element of the four fermion operator to a 
product of bilinear matrix elements. The bilinears are either the time component of the axial current 
or the pseudoscalar density. Each of these bilinear matrix elements is separately averaged over the 
gauge configurations. In the following we generically refer to the ratios of matrix elements as B- 
parameters. Using ratios cancels both the exponential decay factors and the overlap of the source 
operators with the kaons. This leads to simplified fitting, and to a reduction in statistical errors. 
For the renormalized AS = 2 LL operator, the I?-parameter one obtains (using the axial current 
in the denominator) is proportional to Bk- To calculate Bj and B 8 and B$ is more complicated, 
since, as explained below, the denominator involves both axial and pseudoscalar bilinears. 

In several applications we need to consider the matrix element of the four fermion operator at 
intermediate stages of the calculations. This we obtain by multiplying the ratio by the product of 
bilinear matrix elements, themselves calculated from 2-point functions. 

A possible source of systematic error is contamination from excited states. In order to assess 
the size of this error, we calculate the ratios with sources having two Lorentz structures: either 
all are pseudoscalars, J ~ P = (ipi^^), or all are axial vectors J ~ A4 = (^17475^2)- The 
resulting 5-parameters should be same in both cases, but we find that the convergence is from 
opposite directions, as shown in Figs. ^| and |3|. If the two asymptotic results overlap then the effect 
of contamination is smaller than the statistical errors. We find that the two results do overlap, as 
illustrated by the Figures, in all channels. The difference between the two cases increases as the 
quark mass is decreased, growing to ~ lcr for the lightest combination, U3U3. In practice, we take 
the average of the results from the two sources as our best estimate. 
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K (p=o) K (p) 




Figure 1: The ratio of correlations functions that we calculate. The sources J are the same for all 
four pseudoscalar mesons. We show the picture for the case of propagation forward in time. 
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Figure 2: (using one-loop matching with /x = g* = 1/a) as a function of the time slice T op at 
which the operator is inserted for the two choices of source J. The data are for the mass combination 
U2U3 with momentum transfer p = (0, 0, 0). The fits are shown by the solid line and the errors by 
the dashed lines. The fit range, T op = 11 — 21, is chosen to be midway between the two sources. 
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Figure 3: Same as Fig. ^| except that the data are for momentum transfer p = (2,0,0). Because 
the signal in the backward propagating correlator (which has the non-zero momentum) dies off for 
T < 15, the fit range is asymmetric. 
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3.3 Statistical Errors 



We estimate errors using single elimination jackknife. We first average the data from the forward 
and backward propagation of kaons on each configuration. Then, within the jackknife loop, we fit 
to the ratios for each of the two sources, and average the results of these fits. 



4 Operator Mixing and Perturbative Running 

In phenomenological applications the matrix elements of operators must be combined with Wilson 
coefficients. These coefficients have been calculated in continuum perturbation theory in the NDR 
scheme. Results are available at a number of scales, and we take 2 GeV as our standard. We thus 
need to relate our results for lattice matrix elements to those for continuum operators in the NDR 
scheme at 2 GeV. 

The general form of the matching formula is 

0^(ji) i = Z(»,a) ij O v * t (a) j , (4.1) 

where % and j label operators. Z can be expanded as a power series in g 2 , and typically only the 
0(g 2 ) term is known. Rather than use the direct perturbative expansion, however, it is better 
to use a renormalization group improved expression. In Appendix [A] we give such an expression 
for Z [Eq. (|A.8|) ] by combining the exact perturbative result of Ji |l8j with the improved lattice 
perturbation theory of Lepage and Mackenzie [HJ. The additional input needed is the continuum 



two-loop anomalous dimension matrix. The improvement is significant numerically when the dif- 
ference between continuum scale ji = 2 GeV and the lattice scales \ j a — i\ j a becomes substantial, 
and when the anomalous dimensions are large. The improvement is most important here in the 
calculations of Bj^. 

In practice we use a slightly different improved expression, that from "horizontal matching", 
which is given in Eq. ( |A.14[ ). The numerical difference between the two methods is smaller than 
our statistical errors. 

We need the matching coefficients for the operators with Dirac structure 7 M L <g> 7 M L (for Br), 
7 M L ® 7^/? (for B 7 and B 8 ), and L ® L (for B s ). Those for the first two structures have been 



previously calculated in Refs. fl20| , |21|, using the DRED (dimensional reduction) scheme or a 
variant thereof. We have extended these calculations by (i) repeating them for all Lorentz scalar 
four-fermion operators, including those needed for Bs, and (ii) matching to NDR in the continuum 
rather than DRED. We give the results in Appendix |B|. We find that using tadpole improved lattice 
operators substantially reduces the one-loop corrections [19L p2| . For the local operators we use, 



tadpole improvement involves changing the normalization of quark fields by replacing y2n with 
y/8n c ^Jl — 3k/4k c . We make this replacement in our numerical simulations, even though it cancels 
between the numerator and denominator of the B parameters. 

It turns out that the results for the matching coefficients for local four-fermion operator can be 
easily generalized from Wilson fermions to those for any improved fermion action. One only needs 
to know the matching coefficients for the five bilinear operators. We present the results in Appendix 
B in such a way that this generalization can be straightforwardly implemented. 

To carry out the matching we need to choose not only the final continuum scale /i, but also a 
scale q* related to the truncation of perturbation theory (and explained in detail in Appendix |A|). 
For our final results we take /i = 2 GeV and vary q* in the range 2 GeV < q* < rc/a, taking 1/a 
for our best estimates. At intermediate stages it is convenient to take fi = q* = 1/a, which we call 



S 



the TAD1 scheme. This corresponds to doing horizontal matching at 1/a, but no further running 
in the continuum. This choice has the advantage that our results depend very weakly on the lattice 
spacing in physical units — the only dependence arises when we have to interpolate to the physical 
kaon mass. 

Finally, we explain how we determine the MS coupling constant, which is needed in the matching 
formula. We do this following Ref. by first solving 



4-7T 

log □ = — oy(3.41/a)(l - 1.185ay) (4.2) 
3 



for ay, and then converting to the MS scheme using 

0555(3.41/0) = a y (e 5/6 3.41/a)(l + 2a v /n) . (4.3) 

For the plaquette we use □ = 0.5937. We run the coupling to other scales using the two-loop 
/^-function. This results in a s {q*) = 0.2049, 0.1927, 0.1523 and 0.1343 for q* = 2 GeV, 1/a, 2/a 
and ir /a, respectively. 



5 Results for Bk 

The kaon S-parameter is defined in the continuum by 

Br{ s = (k \q\k- ) _ ^ 

(8/3) (K \s a ^ v j 5 d a \0) (0\s b j v ^ 5 d b \K ) ' 

where Q is the left-left AS = 2 operator 

Q = Sa'j^Lda Sb'j^Ldb , (5.2) 

and L = (1— 75). We are ultimately interested in the value of Bk when Q does not insert momentum 
and when the "kaons" have physical masses. The definition Eq. [TT] is, however, more general, and 
is useful at intermediate stages. B K depends on the renormalization scale, /1, because Q has a 
non- vanishing anomalous dimension unlike the axial current. 

To calculate Bk, we use the results of Appendices [A] and [FJ to match Q to a linear combination 
of lattice operators with different Lorentz structures. Thus Bk itself is a linear combination of 
lattice S-parameters, which we define to be the ratio of the matrix element of the four-fermion 
operators [listed in Eq. flB.2|) 1 to 

(4/3)(ir |(s a 7,75^) lat |0)(0|(s fe7i ,7 5 4) lat |^o). (5.3) 
The four-fermion operators in Eq. ( p.2|) have half the Wick contractions of Q, which is why Eq. ( p.3[) 



contains 4/3 in contrast to the 8/3 in Eq. (|5.1| ). Note that we do not include perturbative matching 
factors in these 5-parameters, i.e. the operators in the numerator and denominator are bare lattice 
operators. In Table [l] we give our results for all the 5-parameters for the U2U3 meson. This is the 
lattice meson whose mass lies close to that of the physical kaon. 

We also include in Table [l] results for Bk using one-loop matching in the TAD1 scheme, i.e. 
fi = q* = 1/a. At tree level, Bk is just the following sum of the ^-parameters (for definitions see 
Appendix 0) 

BK-Al + Al + Vl + Vl + Al + A^ + Vl + Vl (5.4) 
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0.28(6) 


B 7 


0.60(0) 


0.60(1) 


0.61(1) 


0.60(1) 


0.60(1) 


0.61(2) 


0.60(2) 


0.62(2) 


0.64(2) 


0.62(3) 


B 8 


0.82(0) 


0.82(1) 


0.82(1) 


0.82(1) 


0.82(2) 


0.82(2) 


0.81(3) 


0.84(3) 


0.86(3) 


0.84(4) 


Bt 


0.59(0) 


0.59(1) 


0.60(1) 


0.59(1) 


0.61(2) 


0.61(2) 


0.61(2) 


0.61(2) 


0.61(3) 


0.60(3) 


Bt 


0.78(0) 


0.78(1) 


0.79(1) 


0.78(2) 


0.81(2) 


0.80(3) 


0.80(3) 


0.82(3) 


0.81(4) 


0.80(4) 



Table 1: i5-parameters (defined in the text) for individual operators without any renormalization 
factors. The results for B^, Bj^ 2 , B^ 2 , and B§ are in the TAD1 scheme. All data are for the 
mass combination U2U3. The data are shown separately for the two types of source J and for the 
five momentum transfers. The columns are labeled by J\p 2 ]J. 



10 





(0,0,0) 


(1,0,0) 


(1,1,0) 


/ -1 -1 1 \ 

(1,1,1) 


(2,0,0) 


subtracted 


cc 


0.920(10) 


0.921(11) 


0.924(14) 


0.925(15) 


0.927(16) 




cs 


0.834(09) 


0.838(11) 


0.841(14) 


0.844(16) 


0.847(19) 




CUi 


0.808(09) 


0.811(12) 


0.813(16) 


0.815(19) 


0.817(22) 




cu 2 


0.796(11) 


0.799(14) 


0.797(20) 


0.799(23) 


0.799(27) 




cu 3 


0.792(13) 


0.795(17) 


0.789(25) 


0.792(29) 


0.787(36) 




ss 


0.558(04) 


0.586(07) 


0.605(10) 


0.623(13) 


0.638(16) 


0.746(216) 




0.444(04) 


0.488(07) 


0.518(11) 


0.545(15) 


0.567(19) 


0.730(165) 


su 2 


0.382(05) 


0.437(08) 


0.473(12) 


0.506(17) 


r\ v s~\ c\ f c\c\ \ 

0.532(22) 


0.716(146) 


su 3 


0.338(05) 


0.402(09) 


0.442(15) 


0.478(20) 


0.508(26) 


0.704(135) 




0.244(04) 


0.332(08) 


0.387(13) 


0.433(17) 


0.466(23) 


0.716(115) 


u ± u 2 


0.112(05) 


0.237(09) 


0.310(15) 


0.369(21) 


0.408(28) 


0.703(098) 


T T T T 

E/1E/3 


r\ r\r\r\ / r\T\ 

0.000(07) 


0.161(10) 


0.248(18) 


0.318(25) 


0.363(34) 


0.687(088) 


u 2 u 2 


-0.091(07) 


0.103(11) 


0.204(18) 


0.283(26) 


0.330(38) 


0.688(082) 


u 2 u 3 


-0.287(10) 


-0.017(14) 


0.112(22) 


0.212(33) 


0.256(49) 


0.670(075) 


u 3 u 3 


-0.597(16) 


-0.191(21) 


-0.021(31) 


0.112(45) 


0.140(68) 


0.643(074) 



Table 2: Lattice results for the i?-parameter of the AS = 2 operator Q in the TAD1 scheme. The 
momentum inserted by the operator is indicated at the top of the columns, and the rows label the 
meson's quark combination. Results after momentum subtraction (see text) are given for the ten 
lightest mesons. 

At one-loop, however, all Lorentz structures enter. In the result for Bk we include one-loop matching 
factors not only in the numerator but also in the denominator, i.e. the lattice axial currents in the 
denominator have each been multiplied by Za = 1 + X(q*)cA- To reduce errors, we always form the 
appropriate linear combination of operators in the numerator and denominator before fitting. 

The results in the Table indicate the relative size of the matrix elements of the contributing 
operators. The largest are the those of the pseudoscalars and tensors. The results using the two 
sources ( J = A4 and P) agree within errors in most cases. Our best estimate, as explained above, 
is obtained by averaging the two results. The errors increase as the momentum inserted by the 
operator increases, but we retain good statistical control for all five momenta. This allows us to see 
clearly that most of the 5-parameters have substantial momentum dependence. 

In Table || we give our results for E>k for all mass combinations. To extract the physical value 
of Bk we must account for the breaking of chiral symmetry by Wilson fermions. The general form 
of the chiral expansion of the matrix elements of four-fermion operators is 

W(j>j) O\K ( Pi )) 



V3)/ 2 



a + (3m K + >ypi ■ p f + 5im K + 5 2 m K pi- p f + 5 3 (pi- p f ) +.... (5.5) 



K,phys 



where /i^phys is the physical kaon decay constant (in lattice units). We are ignoring chiral loga- 
rithms and terms proportional to (m s — m^) 2 . The former are difficult to distinguish numerically 
from the terms we include, while the latter we expect to be small, especially for the range of quark 
masses studied. We are also assuming that we do not need to include possible p 4 terms which 
violate Lorentz symmetry. 

For most operators we expect all terms in Eq. ( |5.5| ) to contribute. For O = Q cont ; however, 
chiral perturbation theory predicts that the matrix element vanishes when Pi ■ Pf = 0, requiring a, 



11 



13 and 8\ to be zero. This means that for Bk, one expects 



fe -SS-¥ ,+ ^ +to " ) ' (5 ' 6) 

With Wilson fermions, however, chiral symmetry is broken explicitly. This would not be a problem 
if one could match with the continuum operator using a nonperturbative method, and if one could 
extrapolate to the continuum limit. We can do neither of these, and so expect non-zero values of 
a, j3 and 8\ proportional to ag 2 and to g 4 . There will also be artifacts of this order contained in 
the other coefficients in the chiral expansion. 

To show that this is indeed a problem, we plot our results for (K°\ Q cont |-^°)/[8/j^ phys /3] in 
Fig. (|. The figure also shows the result of fitting the data with the six-parameter form of Eq. (|5.5| ). 
A good fit is obtained, with non-zero values for the artifacts a and (3. The curvature in the data 
indicates the presence of p A terms, but the values of the individual <5« are poorly determined. 

The curvature is drastically reduced if we consider the ratio 



(K°\Q\K<>) 



mi 



Bjcpi-pf, (5.7) 



as shown in Fig. |5|. Here fx is the lattice result for the decay constant for the corresponding meson 
|22|| . This simplification has been noted previously by the APE group J7|. The figure also shows 
a fit to the same form as on the r.h.s. of Eq. |5.5| . Not only are the higher order terms smaller 
(although they remain poorly determined), but also (3 is much reduced. Fitting to this ratio leads 
to smaller errors in the final result, and we use it in all of the following analysis. 

We can reduce the uncertainty in the higher order terms by noting that the matrix element of 
Q cont /[8/|-/3] is to good approximation a function only of pi ■ pf, as shown in Fig. |6]. This feature 
was also found in Ref. @. This suggests fitting with 5i = 62 = 0, and the result of such a fit is 
shown in the figure. Thus the dominant artifact is a, with a small contribution from j3. that leads 
to the observed deviations from a single curve. This fit leads to the smallest errors in B^, and we 
take it as our standard. 

In order to extract a value for Bk we use a procedure similar to that we advocated in Ref. [|J, 
except that we now have better control over the higher order terms. Having obtained a fit to the 
data, we simply discard the terms which we know to be artifacts. Applying this to the fit shown in 
Fig. H gives 

B K (fx = q* = l/a) « 7 + {5 2 + 6 3 )M K>phys = 0.68(4) . (5.8) 

We obtain consistent results, with larger errors, from the fits of Fig. | [0.68(6)] and Fig. | [0.65(10)]. 
The subtraction procedure has, however, increased the errors substantially compared to those in 
the values for individual mass combinations. For example, the error for U2U3 in Bk is 0.01 (see 
Table g). 

To give an idea of the importance of the p A terms we have done a variety of alternate fits, all 
to the ratio Q/[8f K /3]. Most interesting are those using only the three lowest order terms (a, (3 
and 7). Fitting to the results from the lightest four mesons alone gives Bk = 0.73(4), while fitting 
to the lightest six gives 0.74(4). These are higher than our preferred number, but the difference 
is only slightly more than la. By contrast, if we fit to the combinations UiU%, SUi and SS (the 
heaviest five of the mesons in Figs. |] and ||), we find a result, 0.80(4), significantly higher than our 
preferred number. Previous work by other groups with Wilson fermions J7L |5|], including our own [|J, 
used three-parameter chiral fits to mesons in this mass range. Our results indicate that the neglect 
of higher order terms might lead to 10-20% systematic errors in their results. 
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a = 


-0.04£ 


3(01' 




^ = 


0.046 


(44) 




- 7 = 


0.668 


(44) 




5 i = 


0.17( 


1.33) 


1 




0.1 1 ( 


2.38 




- 6 3 = 


0.51( 


1.07) 




m 



K 



Figure 5: Results for matrix element Q cont / [8/^/3], together with the results of a chiral fit. 
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.1 



a =-0.046(001) 



7 = 0.657(046 
6 ? = 0.526(077 





o p 

x p 

x p 

x p 

X p 



111 

110 
000 



0.1 0.2 

Figure 6: Results for matrix element of Q cont / [8/|-/3], together with the results of a chiral 
excluding Si and S 2 - The data are plotted as a function of Pi ■ pf. 
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We have checked our result by repeating the analysis using a second method. Returning to the 
original definition of the chiral expansion, Eq. ( |5.5| ), we note that contributions from the artifacts 
a, (3, 5i can be canceled by combining pairs of points at different momentum transfers: 



Here E^ is the energy of the lattice kaon with momentum The r.h.s. of Eq. ( |5.9|) differs from 
the form for Bk, Eq. ( |5.8| ), by 5%M(Ei + E 2 — M). We use the value of £3 extracted from the fit 
in Fig. |5| to correct for this. The results of this analysis for the ten light mass combinations are 
given in the sixth column of Table |2|. We have averaged the results from [g 1; q 2 ] — [(100), (000)] and 
[<?i)<?2] = [(HO), (000)], the channels with the best signal. Extrapolating/interpolating to M K , we 
find Bk = 0.66(9), in good agreement with our earlier result. 

We thus are confident that our results do not have significant errors due to yet higher order 
terms (of 0(p 6 )) in the chiral expansion. What is more difficult to gauge, however, is the size of the 
lattice artifacts in the chiral coefficients which we keep, i.e. 7, 5 2 and £3. These artifacts are due 
both to discretization errors and to higher order terms in the perturbative matching coefficients. 
Recently the APE group has studied Bk using Wilson fermions and the improved "clover" action, 
which removes errors of 0(ag 2 ) 0. Their result suggests that the largest source of error in our 
calculation is from the matching coefficients rather than discretization. 

Assuming that truncated perturbation theory is the problem, we can estimate the size of the 
associated systematic errors using a method adapted from that of Bernard and Soni M. As explained 
in Appendix |B|, Q matches with a linear combination of five lattice operators: itself (20f) and the 
four other Fierz self-conjugate operators 2 _ 5 . The idea is to separately adjust the matching 
coefficients of 2 _ 5 relative to that of Of so as to remove the artifacts a, (3 and 8±. The change 
in matching coefficients would lead to a change in the value of Bk- Since there are four operators 
to adjust but only three coefficients to fix at 0{j) A ), we need to make additional assumptions to 
carry out this procedure. One approach might be to consider three operators at a time, and use 
the spread in the final results as an indicator of the systematic error. 

In practice, we have only been able to carry out this procedure semi-quantitatively because our 
statistical errors are too large. In table |3| we show the results of chiral fits to matrix elements of 
various operators keeping only a, 7, and 5 3 : 

1. The perturbatively corrected AS = 2 operator Q (in the TAD1 scheme) divided by 8f K /3. 
This is the matrix element relevant for Bk, the results for which are shown in Fig. |^. 

2. Four linear combinations of the bare lattice operators 2 _ 5 , each divided by 8/^-/3. 

The three parameter fit requires that the data lie on a single curve when plotted against Pi ■ pf, 
which is not a bad representation of the Q data, as shown in Fig. |6]. The fit to 2 _ 5 is better still. 
We do not introduce the fourth parameter f3, because, although this improves the fit for Q, it makes 
little difference for 2 _ 5 , and increases the errors in the fit parameters. 

Since our data is reasonably represented with only the single artifact a, and since the four 
"subtraction operators" can be similarly represented, our subtraction procedure is as follows. Take 
each of the four operators listed in the table in turn, add them to Q so as to cancel a, and use the 
change in J+S3m KphyB as the estimate of the possible shift in B K - In fact we ignore the contribution 
of £3, since it is both small and uncertain. The shifts in Bk are then 0.03, 0.07, 0.06, and 0.05 for 
the four subtraction operators listed in the Table. Thus we take 0.05 ± 0.05 for our estimate of this 
shift, the generous error accounting for the uncertainty in this subtraction method. Thus our final 
estimate of B K is 



E 1 B K (qi) - E 2 B K (q 2 ) 
(£1 - E 2 ) 



7 + 5 2 M 2 + 5 3 M(E 1 + E 2 ). 



(5.9) 



B K ((j, = q* = 1 1 a) = 0.73 ± 0.04(stat) ± 0.05(subt) . 



(5.10) 
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Operator 


a 


7 


Ss 


Q 


-0.0457(14) 


0.713(20) 


0.52(7) 


ot + ot 


1.79(5) 


1.1(5) 


-1.5(1.4) 


ot - ot 


-3.36(9) 


-4.8(9) 


2.9(2.9) 


2(Of + Ot) 


-1.29(4) 


-1.8(4) 


0.3(1.2) 


2(Of - Ot)/3 


-0.71(2) 


-0.8(2) 


0.0(6) 



Table 3: Results (in lattice units) of three parameter chiral fits to various operators. 

Running in the continuum to fi = 2 GeV increases the result slightly 

B K (n = 2 GeV, q* = l/a) = 0.74 ± 0.04(stat) ± 0.05(subt) . (5.11) 

We close this section with some comments on the reliability of our result. 

1. Our result is almost independent of the horizontal matching scale q*: it increases only by 0.02 
if we use q* = ix/a. This is in apparent contradiction with the conclusions of Ref. 0, who 
find a large dependence of Bx on their choice of a s . The difference is due to the fact that we 
have used tadpole improved perturbation theory, in which the perturbative corrections to both 
numerator and denominator of Bk are small. Although the tadpole improvements formally 
cancel at one-loop order, there is a significant numerical difference due to 0(a 2 ) terms. 

2. The changes to the coefficients of the subtraction operators required to cancel the a term in 
the chiral expansion are smaller than the one-loop perturbative results. For example, we can 
cancel a by increasing the coefficient of Ot + Ot from — 1.0a s /7r (its one-loop perturbative 
value) to — 0.6a s /7i ~ —a s j-n + 7(« s /7r) 2 . This is not an unreasonable value for a higher order 
perturbative term. 

3. Ideally, one would like to determine the matching coefficients non-perturbatively. First results 
with Clover action |7|, |8[], and with Wilson fermions @ have been obtained. The latter are 
consistent with our result obtained using the chiral fits and perturbative matching. 

4. Finally, we note that the matrix elements of V 2 and V 1 , although they are the largest of all 
Lorentz structures (see Table |T]), make only a small contribution to Bk- This is because they 
depend, to good approximation, only on m 2 K and not on p^ ■ pf as exemplified in Fig. [F| On 

the other hand, for B^ 2 and B^ 2 , where V 2 and V 1 are part of the continuum operator, these 
operators give the dominant contribution. 



6 B D 

The analysis of our results for 

B _ (Dolcg-ff.Lug c b ^^Lu h \Do) 

(8/3) (A)|e a 7*/y5«a|0) (0 1 Q,7v75«6 1 D ) ' 

and Bd s (related by u — > s) is more straightforward. There are no significant constraints from 
chiral symmetry — all terms in the expansion of Eq. ( p.5| ) are expected to be present in the physical 
matrix element. Thus we simply calculate the S-parameters and extrapolate or interpolate to the 
physical light quark masses. 
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Figure 7: Three parameter fit to the matrix elements of V 2 / [8/|-/3]. The data for V 1 is similar and 
the VSA relation V 2 = 3V 1 is approximately valid. 
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Our results for CS and CUi mesons are given in Table [|. They show no significant variation 
with momentum transfer, and very little dependence on the light quark mass. Thus the a term in 
Eq. ( |5.5D is dominant — in striking contrast to the results for light-light mesons. Quenched chiral 
perturbation theory predicts that Bo should diverge as the light quark mass vanishes p3| , [24]], with 
the effect beginning at m SiP h ys /3 — m SiP hy S /4. We see no sign of such behavior with our light quark 
masses, which extend down to approximately m SiPhys /3. 

To investigate the convergence of perturbation theory, we show, in Table |4], separately the "di- 
agonal" and "mixing" contributions to Bp. The former is the contribution from the lattice operator 
having the same form as the continuum operator (i.e. 20±), while the latter is the contribution from 
the other operators introduced by one-loop mixing (C^s). The mixing contributions are ~ 15% of 
Bp, indicating reasonable convergence. Their fractional contribution is also almost independent of 
the momentum and the mass of the light quark. 

Our final results are 

B D (fi = q* = l/a) = 0.777(15) and B Da (n = q* = l/a) = 0.815(9) . (6.2) 
These convert to 

B D {p, = 2 GeV, q* = l/a) = 0.785(15) and B D ,(/j, = 2 GeV, q* = l/a) = 0.822(10) . (6.3) 

The dependence on q* is, as for Bk, smaller than the statistical errors. Since Bo and Bo s are 
correlated, their ratio is well determined 

B Ds /B D = 1.048(15). (6.4) 

The variation of Bo with the light quark mass is described by 0.777(15) + 0.77(27)m 9 a. Here m q a 
is the quark mass obtained from the Ward identity for the axial current, and m s is fixed using 
ill. 



7 Results for Bj /2 and B 3 8 /2 

We now turn to the left-right operators which appear in the effective weak Hamiltonian due to 
electromagnetic penguin diagrams: 

Q-7 = {saj^Lda) [(ubjftRub) - \{d b ^^Rd b ) - \(s b ^^Rs b )} , (7.1) 
Qs = (salnLd b ) [(ublpRua) - l{db^Rd a ) - \{'§ b ^ ll Rs a )} . (7.2) 

We have results only for the 1 = 3/2 parts of these operators, since the 1 = 1/2 parts have penguin 
contractions which we have not calculated. Fortunately, the 1 = 3/2 parts are of phenomenological 
interest, because they are the only operators giving rise to an imaginary part in the K + — > 7r + 7r° 
amplitude, and thus give the dominant electromagnetic contribution to e'/e. We calculate B- 
parameters for K + — > tt + matrix elements, which are simpler to calculate than those for K — > tttt, 
but are equal at leading order in chiral perturbation theory. Thus our 5-parameters are defined in 
the continuum by 

7 l(7i+\u 75 d\0)(0\s 75 u\K+) - (7c+\u lm u l5 d\0)(0\s 1 ^ 5 u\K+) ' 1 ' ; 
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P= (0,0,0) 


p = (i,o,o; 


V 


= (1,1,0) 


P= (1,1,1) 


p=(2,0,0) 


cs 


Full 


0.83(1) 


0.84(1) 




0.84(1) 


0.84(2) 


0.85(2) 


Diagonal 


0.71(1) 


0.71(1) 




0.72(1) 


0.72(1) 


0.72(2) 


Mixing 


0.12(0) 


0.12(0) 




0.12(0) 


0.12(0) 


0.12(0) 




Full 


0.81(1) 


0.81(1) 




0.81(1) 


0.82(2) 


0.82(2) 


Diagonal 


0.68(1) 


0.69(1) 




0.69(1) 


0.69(2) 


0.69(2) 


Mixing 


0.13(0) 


0.13(0) 




0.13(0) 


0.12(0) 


0.12(0) 


cu 2 


Full 


0.80(1) 


0.80(1) 




0.80(2) 


0.80(2) 


0.80(3) 


Diagonal 


0.67(1) 


0.67(1) 




0.68(1) 


0.68(2) 


0.68(2) 


Mixing 


0.13(0) 


0.13(0) 




0.13(0) 


0.12(0) 


0.12(0) 


cu 3 


Full 


0.79(1) 


0.79(2) 




0.80(2) 


0.79(3) 


0.79(3) 


Diagonal 


0.67(1) 


0.67(2) 




0.67(2) 


0.67(3) 


0.67(3) 


Mixing 


0.13(0) 


0.13(0) 




0.13(0) 


0.12(1) 


0.12(1) 



Table 4: Results for B D (in the TAD1 scheme) for the heavy-light mesons CS, CUi as a function of 
momentum insertion. Also shown are the separate contributions of the diagonal and mixing parts 
of the operator. An error value of means that it is less than 1/2 in the last significant digit. 
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TTl 11 

Full 


B 7 
Diagonal 


mixing 


TTl 11 

Full 


^8 

Diagonal 


mixing 


cc 


0.859(09) 


1.07(1) 


-0.22(0) 


0.941(10) 


1.23(1) 


-0.30(0) 


cs 


0.822(09) 


1.08(1) 


-0.26(0) 


0.943(10) 


1.22(2) 


-0.29(0) 


CUx 


0.812(10) 


1.07(1) 


-0.27(0) 


0.940(11) 


1.20(2) 


-0.28(0) 


cu 2 


0.806(11) 


1.06(2) 


-0.27(0) 


0.936(12) 


1.20(2) 


-0.28(1) 


cu 3 


0.800(13) 


1.06(2) 


-0.27(1) 


0.932(14) 


1.19(3) 


-0.28(1) 


ss 


0.723(05) 


1.05(1) 


-0.34(0) 


0.917(06) 


1.17(1) 


-0.26(0) 


SUx 


0.695(04) 


1.02(1) 


-0.34(0) 


0.900(05) 


1.14(2) 


-0.25(0) 


su 2 


0.680(04) 


1.01(2) 


-0.35(1) 


0.890(06) 


1.12(2) 


-0.25(0) 




0.669(05) 


1.00(2) 


-0.34(1) 


0.880(06) 


1.11(2) 


-0.25(1) 




0.660(04) 


0.99(2) 


-0.35(1) 


0.874(05) 


1.10(2) 


—0.24(1) 


u x u 2 


0.640(04) 


0.97(2) 


-0.35(1) 


0.858(05) 


1.08(3) 


-0.24(1) 


T T T T 

UIU3 


0.625(04) 


0.96(3) 


—0.34(1) 


0.843(05) 


1.06(3) 


r\ no/1 \ 

—0.23(1) 


u 2 u 2 


0.617(04) 


0.95(3) 


-0.34(1) 


0.837(05) 


1.05(3) 


-0.23(1) 


u 2 u 3 


0.599(04) 


0.93(3) 


-0.34(1) 


0.819(06) 


1.03(3) 


-0.23(1) 


u 3 u 3 


0.578(05) 


0.91(4) 


-0.34(1) 


0.797(07) 


1.01(4) 


-0.22(1) 



Table 5: Results for B 7 and B g in the TAD1 scheme (/i = q* = 1/a). 

The operators Q 7 ^ are defined in Eqs. ( |B.25| ) and ( |B.26| ). 

The operators in are matched onto lattice operators using the results of Appendix ||. 
A slight complication arises because the denominators consist of two terms, one of which (the 
pseudoscalar) depends on the renormalization scale /i, while the other (axial) part does not. Thus 
we cannot simply use the 5-parameters defined above which have only the axial-currents in the 
denominator. Instead, we use horizontal matching followed by 2-loop running to calculate both 
numerator and denominator at the final scale /i, and then take the ratio, all within the jackknife 
loop. We repeat the entire procedure for various choices of q* and \i. 

The calculation involves no chiral subtractions, since all terms in the chiral expansion are phys- 
ical. To show the quality of the results we display B 7 and B^ 2 (in the TAD1 scheme), in Figs. || 
and H respectively, for the U 2 U 3 meson. The statistical errors are small, and the results using the 
two sources agree. 

We have done the calculation of these matrix elements only for the case in which the initial 
and final mesons have the same quark composition. Our results in the TAD1 scheme are collected 
in Table [5[ We see the same gradual increase towards unity as the meson mass increases that we 
observed between B^ and Bp. Diagonal and mixing contributions are defined as for Bp. We have 
also included results at different momenta for U 2 U 3 mesons in the last two rows of Table 0. The 
lack of momentum dependence is due to the dominance of the pseudoscalar matrix elements, which, 
as noted above, depend only on m 2 K . 

In Table || we give results for fi = 2 GeV with various values of q*. There is a significant 
variation with q*, considerably larger than the statistical errors. This indicates that the perturbative 
expansion for the anomalous dimensions is converging more slowly than for Bk and Bo- We take 
q* — 1/a for our final result, and use the range of values for q* = 2 GeV — n/a to estimate the 
systematic error. This leads to 

B 7 /2 (NDR, 2 GeV) = 0.58 ± 0.02(stat) +0 '° 7 (pert) , (7.5) 
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Figure 8: Results for the ratio of correlators defining B 7 ' , in TAD1 scheme for the U2U3 meson 
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Figure 9: Results for B^ 2 in TAD1 scheme for the U 2 U 3 meson. 
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2 GeV 0.549(18) 0.790(27) 

1/a 0.578(19) 0.807(27) 

2/a 0.641(21) 0.835(28) 

Tt/a 0.654(21) 0.837(28) 



Table 6: Variation of B>? and Bg, run to /x = 2 GeV, for various choices of q*. 



Bg (NDR, 2 GeV) = 0.81 ± 0.03(stat) ' (pert) . (7.6) 

Note that the statistical errors are somewhat larger after running to fi = 2 GeV (compare results 
in Table |6] to those in Table |5|). This is because, in order to evaluate the B— parameters at 2 GeV, 
one has to combine results of three separate fits after the running. Consequently, the cancellation 
of errors between the numerator and denominator is not as good. 



8 Results for Bs 

We close with results for the ^-parameter introduced in Ref. to calculate the difference of 

B-meson lifetimes _ _ 

B =B + = (Bs\bgLs a b c Ls c \B s ) = 

S ~ 4 (5/3)(B s \b al5 s a \0)(0\b al5 s a \B~) ' 
We also consider the operator related by a Fierz transformation 

B+ (B s \b a Ls c b c Ls a \B s ) 

5 ~ (-l/3)(B s \b al5 s a \0)(0\b al5 s a \B;) ' 1 ' j 

The perturbative matching coefficients for this operator have not been calculated previously — they 
are given in Appendix |||. The 1-loop improved operators Of and Of are defined in Eqs. B.36j and 
B.37| . (One needs the matrix elements of both operators due to the mixing between them.) The 



calculation of their matrix elements and the corresponding B-parameters is straightforward, and 
the results for the 15 mass combinations, as a function of the matching scale /i = q*, are given in 
Tables [7] and || The results after extrapolation of the light quark to the physical strange quark 
mass are given in Table [| 

Note that, unlike the previous cases of .B-parameters, we cannot run these results to a final 
common scale. As can be seen from the results of Appendix A, this requires knowledge of the 
two-loop anomalous dimensions. It would be inconsistent to use one-loop matching and then only 
run using one-loop anomalous dimensions. Thus we cannot study the q* dependence of the results 
at fixed final scale \i. The fact that both £> 4 + and B 5 +, evaluated at \i = q*, depend only weakly 
on q*, as is apparent from the Tables, is an accident. The origin of this feature can be traced to 
the large value of the matrix elements of the lattice mixing operators, and consequently significant 
variation with a s (q*). 

The best we can do at present is to use our preferred value of q* = 1/a, and quote the results 
only for \x = q*, i.e. using horizontal matching. These results are 

5+(NDR,/i = q* = 1/a) = 0.80 ± O.Ol(stat) , (8.3) 
S+(NDR,ii = q* = 1/a) = 0.94 ± O.Ol(stat) . (8.4) 
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Table 7: The B-parameter for the operator Oj . The data are shown for all 15 mass combinations, 
and for various choices of ji — q*. 
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Table 8: The B-parameter for the operator O5". The data are shown for all 15 mass combinations, 
and for various choices of ji — q*. 
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Table 9: The B-parameter for the operator 0\ and O5 after horizontal matching to the continuum 
at scale q*. The CUt data have been extrapolated to rris^M^) in the light quark mass. 
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A Issues in matching lattice and continuum operators 

In this appendix we discuss the matching of lattice and continuum operators. We also collect results 
for the anomalous dimensions needed for the particular matrix elements we calculate. 

We wish to find the lattice regularized operators which match with (i.e. have the same physical 
matrix elements as) the continuum operators of interest: 

O cont (/i) = Z(ji, a)C latt (a) . (A.l) 

Here O cont is a vector of operators which is closed under continuum and lattice mixing, and Utt 
the corresponding vector of lattice operators. By comparing matrix elements, one can calculate the 
matrix Z in perturbation theory. The one-loop result has the form 

Z(n,a) = l + Xc(fia), (A.2) 

where A = g 2 /(!6ir 2 ) and c is the one-loop matching matrix. Appendix [B| describes the calculation 
of c for four-fermion operators. The issue we address here is the choice of g to use in A. Should it 
be the continuum coupling fl^jg , (/x), or the bare lattice coupling g(a), or some sort of an average? 
The uncertainty this introduces can be comparable to the statistical error, as is the case in the 
calculation of Br using staggered fermions |jTT| . 

To address this issue we use the exact perturbative formula for Z given by Ji |L8]. Generalized 
to the case of mixing, this is 

rr I \ rp ( f 9M , „ 7cont(ff") ^ „ ( f° , /7latt(fl0\ / \ Q\ 

Z(ji, a) = T g n exp - / dg — - T g , exp - / dg - — -— . A.3 

The 7's are the anomalous dimension matrices, which are different in the continuum and on the 
lattice. The /3's are the usual /3-functions, again different in the continuum and on the lattice, 
though only at 0(g 7 ). The T 5 's indicate that the exponential integrals are g-ordered. This formula 
has a simple physical interpretation: first use the renormalization group to run on the lattice to a 
lattice spacing as small as one wants, then match with the continuum at tree level, and finally run 
back to the desired scale in the continuum theory. It is exact except for possible non-perturbative 
corrections. 
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To use Ji's formula we have to specify both the continuum and lattice renormalization schemes. 
In the continuum we take NDR, (i.e. MS plus a particular set of rules for dealing with 75 in n 
dimensions) so that is in the MS scheme. On the lattice, Ji uses the bare lattice coupling 
constant, and takes the scale to be a. The work of Lepage and Mackenzie suggests, however, that 
it is better to reexpress all lattice perturbative expressions in terms of a continuum-like coupling 



19]. Thus we suggest improving Ji's formula by using the same definition of renormalized coupling 
constant on the lattice as in the continuum. In other words, we imagine expressing perturbative 
series on the lattice, e.g. that for the anomalous dimension matrix, in terms of gy^{q*), where 
q*a = K, with K a constant. This is just a change of variables 

0ms(?*> 9(a)) = g{a) + Cl {q*a)g{af + 0{g 5 ) , (A.4) 

with c\ a known constant with known dependence on q*a. The value of q*a should be chosen so as 
to improve the convergence of perturbative series. 

To implement this improvement, we change variables in the integral over g' in Eq. ( |A.3| ), from 
g' to g^g(q* = K/a,g'). The result is 

,7iatt0') f° dg m 7iatt(y(#Ms)) 



^ / /latxyj 1 _ 1 -VMS native WMS// g\ 

s(a) Aattty) 4^(g*=if/a) \dg m (g') / dg'\ (-dg'(a)/d\na) 



7iatt(ff'(ffMs) ) 
QMsii^K/a) ^ MS dg^s/ d In q 



dm T Ky ^!: (a.6) 



m (q*=K/a) Pcont{9Ms) 

Thus we end up with the same form, but with the continuum /5-function in the denominator, a 
different lower limit, and the anomalous dimension matrix expressed in terms of g^. This change 
of variables makes no difference if we work to all orders, but does affect the result when we truncate 
perturbation theory. The hope is that 7i at t converges more quickly as a result of this reorganization. 

In practice we know only the two- loop anomalous dimensions for the operators of interest, so we 
expand out Eq. |A.3| (with the substitution |A.7| ) to second order, following the method of Ref. @ 



Zfa, a) = [l + A(ai) Jcont + OiX 2 )}- 1 W- 1 i^j W[l + A(0-WO + 0(X 2 )} . (A.8) 

Note that the g 2 in A is in the MS scheme wherever it appears. Note also that the a depen- 
dence comes in implicitly through q*. 7^ is the result of diagonalizing the the one- loop anomalous 
dimension matrix, 70, (which is the same in both lattice and continuum schemes) 

7o D = W^W- 1 . (A.9) 

The 2-loop contribution feeds in through 



j =^.-W- 1 MW . (A.10) 



where the matrix M is 
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J cont differs from Ji att because the two-loop anomalous dimension matrices differ. In practice, the 
lattice J's contain more off-diagonal terms than those in the continuum. For example, in the 
continuum the operators Q-? and Qs only mix with each other (and with evanescent operators 
which vanish in 4-dimensions) , while Ji att mixes these operators also with Of 45 , as described in 
Appendix |B[ The structure of Eq. |A.8| means, however, that, working from left to right, we can 
stay in the two dimensional space spanned by (Q 7 , Q 8 ) until we encounter J latt . Thus we need only 
the rectangular piece of Ji att which connects (Q7, Q 8 ) to the full basis. 

Equation ( |A.8|) is (aside from a minor modification discussed below) the result which we propose 
using in place of Eq. ( |A.2|) . It makes clear which scale to evaluate the coupling at wherever it 
appears: the lattice scale in the lattice anomalous dimension, and the continuum scale in the 
continuum anomalous dimension. It is valid even if /ia is substantially different from unity, for 
it sums up the large logarithms of /ia. It can be applied to tadpole improved operators — tadpole 
improvement simply changes the two-loop lattice anomalous dimension, presumably making the 
perturbative series converge more quickly. 

To apply Eq. [A.8| we need to know the J's. In the continuum, the two-loop anomalous dimensions 



are known for most of the operators of interest (they are listed below), and from these we can 
construct J con t- Jiatt can then be determined by comparing the one- loop matching result to the 
general formula, in the following way. We recall that, when deriving Eq. ( |A.2|) , one equates the 
lattice and continuum coupling constants (their difference being a higher order effect). For this to 
be true in the general formula we must set fi = q*. Thus we have 

Z(ji = q*,a) = [1 + A(0 c(q*a) + 0(X 2 )] 

= [1 + A(g*) J cont + 0(A 2 )]- X [1 + A(g*) J latt (g*) + 0(A 2 )] , (A.12) 

and so 

Jiatt (g*) = Jcont + c{q*a) . (A. 13) 

We now see explicitly that Ji a tt depends on q*, and upon whether we tadpole improve or not (which 
changes c —>■ c — tad, where tad is the tadpole contribution). 

To express the result at scale /1, we have in fact used a slight variant of Eq. A. 8 , which we call 
"horizontal matching" : 



Z(fi,a) = [l + A^Jeo.t]- 1 ^- 1 (^y) 2 "° W[l + X(q*)J CO nt][l + X(q*)c(q*a)]. (A.14) 

This differs by corrections of size A(g*) 2 , which turn out to be much smaller than our statistical 
errors. The physical interpretation of horizontal matching is that we first match lattice and contin- 
uum operators at the scale q* (using the same coupling for both), and then run, in the continuum, 
down to the final scale \i. The important practical advantage of both horizontal matching and that 
based on Ji's approach is the improved treatment of terms involving A, when A(/x) is significantly 
different from A(g*). We have, in effect, used the renormalization group to include a subset of 0(A 2 ) 
terms. 

The remaining issue is the choice of q*. That we have such a choice simply reflects the fact that 
we have truncated the perturbative expansions for 7 and (3. Lepage and Mackenzie have suggested a 
scheme for estimating q* which works for finite lattice renormalizations jn|. We have not, however, 
found a way to implement their scheme for divergent operators. Presumably one could use the 
alternative BLM scheme ||25|| , but this involves calculating a subset of diagrams for the two- loop 



matching coefficients (or three- loop anomalous dimensions), which has not been done. Thus we 
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have simply carried out the calculation for a reasonable range of q*: 2 GeV, 1/a, 2/ a and n/a. 
Previous work suggests that, for tadpole improved operators, the optimum q* lies towards the lower 
end of this range For this reason we quote our final value using q* = 1/a. It is important to 



realize, however, that this intuition concerning q* may not apply to matching calculations such as 
those we are using here. Thus we use the results from the entire range of values of q* to estimate a 
systematic error. 

We close this appendix by collecting results needed to carry out the matching just described. 
We quote all results for rif = 0, and give anomalous dimensions in the NDR scheme. The operators 
are listed in Eqs. ([R21HR271) . 

The coefficients (3q and (3\ of the (3— function (fidg/dfi = —f3o9 3 — Pi9 5 ~ ■ ■ ■) are independent of 
the renormalization scheme, and are standard. The anomalous dimensions for Q are 



7o 



7i 



-7. 



(A.15) 



The operators Q^ 2 and Q^/ 2 mix, with anomalous dimensions 



7o 
7i 



71/3 
-225/2 



-99 
•1331/6 



(A.16) 



The result for 70 is standard, while we have extracted 71 from the results of Ref. [^] by considering 
only a subset of the graphs required for the full operators 0,7$. The results for the operators needed 
to calculate B s are given at the end of Appendix M. 



B One-loop matching for four-fermion operators 

In this appendix we give a general formula for the one loop matching coefficients for four-fermion 
operators of the form r <S> V . This formula relates continuum operators defined in the NDR scheme 
to lattice operators which are local, i.e. reside on a single lattice site. All one needs to know are 
the matching coefficients for the five bilinears S, P, V, A and T. The formula applies not only 
for Wilson fermions, but for any improved Wilson fermion action, e.g. the Sheikholeslami-Wohlert 
action. It is not useful for staggered fermions, where many of the operators of interest are not local 
(see, for example, Ref. P7fl ). 

Matching coefficients for a subset of four-fermion operators were calculated long ago by Martinelli 
[2H[ ] and Bernard et al. |?T|. These authors use different regularization schemes in the continuum, 



namely DRED (dimensional reduction) and DREZ (a variant of dimensional reduction). One of our 
purposes here is to convert these results to the NDR scheme, which is the standard scheme used in 
the evaluation of continuum coefficient functions. A second purpose is to extend the calculations to 
operators of the form S <8> S + P <S> P and T <g> T. As discussed in the text, the matrix elements of 
one such operator has recently been shown to be of phenomenological interest. Our final purpose 
is to present all the results in a simple form which is easy to evaluate. 

Our method is adapted from that used for staggered fermions in Ref. G7J. It proceeds in two 



stages. The first stage uses the method of Martinelli p0[ , who pointed out that all the one- loop 
vertex diagrams for four-fermion operators can be brought into the form of bilinear corrections 
using Fierz transformations and charge conjugation. This works, however, only if the continuum 
scheme keeps the gamma-matrices in four dimensions, which is why he chose DRED. Following 
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Action 


Scheme 


cs 


cv 


ct 


ca 


cp 


Wilson 


NDR 


-12.952 + 6£ 


-20.618 


-17.01 - 2£ 


-15.796 


-22.596 + 6£ 


Wilson 


DRED - NDR 


1 


0.5 


-1 


0.5 


1 



Table 10: Bilinear matching coefficients, with i = log(/ia). Tadpole improvement (discussed in the 
main text) adds t&d/Cp = 12.86 to each of these numbers. 



Ref. |27 |, we actually use a slightly different intermediate scheme, NDR'. In the second stage we 



convert from NDR' to NDR. This part of the calculation is done in the continuum and it cannot 
be reduced to that for bilinears. It is, however, independent of the lattice fermion action. Much of 
the work for the second stage has been done in Ref. H27f . Our only new result is the construction 
of an n-dimensional definition of the operators S (g) 5* + P <8> P and T <g> T which maintains Fierz 
symmetries. This is not an essential condition, but it is desirable for aesthetic reasons. 

We begin with some notation. The one-loop matching coefficients for bilinears are defined by 

gcont = B Ut (1 + Cp x c . (/m)) ^ ( g ^ 

where i = S,V,T, A, P, Cp = 4/3, A = g 2 /167r 2 , fi is the renormalization scale of the continuum 



operator and a the lattice spacing. The values of the Cj are given in Table |10| for Wilson fermions 
in both DRED and NDR schemes. These values are taken from Refs. |2U|, EH], and converted to 
NDR using the conversion factors given in Ref. [^| . Note that the corrections are uniformly small 



after tadpole improvement For completeness we mention that the analogous results of the 

renormalization constants for bilinears with the Sheikholeslami-Wohlert (clover) action have been 
calculated by Gockeler et al. p9|. 



The four-fermion operators we consider are linear combinations of the five Lorentz scalars 



s = 




V = 




T = 








A = 


^(^17^75^2) (^37/^75^4) , 


V = 


(^175^2) (^375^4) • 



(B.2) 



Keeping the four flavors distinct eliminates "penguin" diagrams, and allows only a single Wick 
contraction for each operator. In the text, we use subscripts to indicate incomplete summation over 
repeated Lorentz indices. For example, A s denotes that the sum over fi is only over 1 — 3 and not 
4, while % indicates that v = 4 in the tensor operator. Each operator also comes with two possible 
sets of color indices, which we distinguish with superscripts, e.g. 

S 1 = (?> 2 b )(^k a ) and 52 = (^>2)(^4) • (B.3) 

The superscript is the number of loops that the color indices form if we take a matrix element of 
the operators between states created by (ip 2 ^' l Pi) an d (i^JPips). 

In the following we use two bases for the Lorentz structures: the "original" basis 

0[ = [S, V, T, A, V] , (B.4) 
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and the "practical" basis (related by a non-orthogonal transformation) 

0l = [(V + A),(V-A), -2(5 -V),(S + V),-(S + V- T)/2] . 



(B.5) 



The advantages of each will become clear in the following. In each basis there is a matrix which 
implements Fierz transformations. This is simplest in the practical basis 



T 



prac 



/ 1 \ 
10 
10 
1 
^00010/ 



(B.6) 



In either basis this satisfies T 2 = 1. There is also a matrix which implements "charge conjugation", 
which here means the effect of transposing the Dirac matrices in one of the bilinears, and conjugating 
the result with C = 7o72- This is simplest in the original basis, 



C OTig = diag( 1,-1,-1, 1, 1) 



(B.7) 



Again, in either basis one has C 2 = 1. Note that this definition is specific to four dimensional Dirac 
matrices. 

To combine the Lorentz and color indices we use the notation 



O-, 



of 



Then, as shown in Ref. pTfl , one- loop matching has the form 



Jat 



OT +A£ (M% C a +M b ij C b +M c ij C c ) of . 



(B.8) 



(B.9) 



There are three types of diagram which contribute, denoted a, b and c, each giving rise to a color 
matrix C and a matrix acting on the Lorentz indices. The color matrices are 
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(B.10) 



The core of the calculation is of the 5x5 correction matrices M. 

The first stage of this calculation is to determine M in a scheme in which the gamma-matrices 
are four- dimensional, say DRED. The "a"-diagrams are just those in which the two bilinears in the 
four-fermion operator are matched independently (see Ref. ||27|| ). Thus M a is simple in the original 
basis: 

M a = 2 diag[c s , c v , c T , c A , c P ] . (B.ll) 



Following Ref. p0| one can then obtain Mb and M c using Fierz and charge conjugation transfor- 
mations: 



Mb 
M c 



TM a T ■, 
-CTM a TC . 



(B.12) 
(B.13) 



These results can then be inserted in the general formula Eq. ( |B.9| ). Fierz symmetry also requires 
that M c be self-conjugate: 

M C = TM C T. (B.14) 
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This result, together with Eq. ( |B.13| ), and the peculiar fact that TCTCT = C, requires that 
[C, A4 a ] = 0. This is trivially true given that M. a and C are diagonal in the original basis. 

This construction fails for NDR, however, because one is then dealing with n-dimensional gamma 
matrices. It is useful, nevertheless, to introduce an intermediate scheme, NDR', which is defined by 
Eq. (|B.11|) , with NDR bilinear matching coefficients, together with Eqs. (|B.12| ) and ( |B.13|) . That 
is, we enforce Fierz and charge-conjugation symmetries by hand. 

The second stage is the extension of the calculation to NDR. For this, we need to define the 
continuation of the four-fermion operators to n-dimensions. There are many ways to do this (cor- 
responding to different choices of "evanescent" operators), and we wish to make a choice which 
maintains the Fierz relations Eqs. ( |B.12| ) and ( |B.14| ). We have not found a way of simultaneously 
maintaining Eq. ( p,13|) , and suspect that this is not possible. Thus, in NDR, Aif, can be obtained 
from M. ai while Ai c is independent. Fierz-symmetric choices for (V ± A) and (S — V) were given 

Definitions for (S + V) and T were given in Ref. 



in Refs. 26 



23, but 



MJ, and we use these. 

the choices made there did not maintain Fierz symmetry. We have found Fierz-symmetric choice, 
which we now explain. 

We work in the "practical" basis. The operators are extended to n = 4 — 2e dimensions such 
that they are the even parity parts of 



0\ = (^7^2) OsTm^) , 

02 = Ol7^2) (^37/^4) , 

5 = -i(l + 3e/4) ffliWuLfo) (faw^ 



(B.15) 
(B.16) 
(B.17) 
(B.18) 
(B.19) 



Here we are implicitly summing repeated indices over all n values, and defining L = (1 — 75) 
and R = (1 + 75). The odd parity parts of these operators are not important for the following 
discussion, but we include them for generality. The analysis is unchanged if we change the sign of 
75. The peculiar factor of (1 + 3e/4) in O5 is necessary in order that, in the final result, hA c is 
Fierz self-conjugate. When we calculate one-loop corrections in NDR the resulting operators have 
the form 

O = diOi + evanescent operators, (B.20) 

i=l,5 

where the evanescent operators vanish when n — 4. We determine the di using projection operators. 
Fierz-symmetric projectors onto 0\^ 3 have been given in Refs. 



26, 30]. To determine d± and d§ 



[up to terms of 0(e )] we use 




l + 5e/12 / (l + e/4) 




(B.21) 



where, in the notation of Ref. [30], the projectors are 



Pa 



— L ® L and 
32 

-2o6^ )27 ^ L07 ^ L - 



(B.22) 
(B.23) 



The matrix in Eq. ( |B.21| ) has been determined by requiring di = 1 if O = Oi, i.e. the basis operators 
project back onto themselves. 
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Using these definitions, we have calculated the matching between NDR and NDR', and thus 
between NDR and lattice operators. We find for our final results 



M a 



( c v+ 

cv- 






Cy-A 
CV+A 










Cs+P 

-\c S -P 



1 



CS-P 



\ 



-2c S -p 

cs+p -1 

\{2c t -c s+ p) 2c T J 



( -( 



V 



c s+ p + 6) 



\cs-p 

jCS-P 



2c S - P 

■{C V+A + Q) Cy^ A 

Cy-A -{cy+A + 6) 

~\{c s+P + 2c T + A) 

-\(c s+ p-2c t + 2) 



2c S -p \ 



\(c s+ p-2ct + 2) 
|(c 5 +p + 2ct + 4) j 



where all the q are in the NDR scheme, and we have defined cs+p = cs + cp, etc. Mb continues 
to be determined from Eq. ( [B.12|) . Note that Ai c is Fierz self-conjugate. The results for NDR' and 
DRED can be obtained by deleting the numerical constants, i.e. keeping only the terms involving 
the Cj. 

We conclude by presenting the explicit results needed for the matrix elements we study in this 
paper. The operators of interest are 



Q 

qT 

Qs 



Ul- 

2 



[(s a 'y ll Ld a )(s b 'j fl Ld b ) + (s a 'y fl Ld b )(s b <y tl Ld a )} , 

{SaJf.Lda) [(UhJ^RUh) - (db^Rdb)] + (SaJ^LUajiubJfj.Rdb) 
(SaJpLdb) [(Uh^^RUa) - (dbJf.Rda)} + (SaJ^LUbjiUblfiRda) 

\[{b a Ls a ){b h Lsb) - |(1 + 3e/4)(6 o 7 At 7 I/ Ls 6 )(& b 7 I/ 7 At Ls )] . 



(B.24) 

^B.25) 
T3.26) 

;b.27) 



Following Ref. [^TJ, we have extended Q and Qs to n-dimensions in such a way that the Fierz 
symmetry is explicit. None of these operators have penguin contractions (in the isospin symmetric 
limit), and so all can be matched onto the lattice using the approach explained above. 

The first step is to rewrite the operators using the four flavor notation. For this purpose it is 
useful to introduce a third operator basis, that of Fierz self-conjugate and anti-conjugate operators: 



Of = [ [o\ ± o\i (o 2 2 ± ol), (ol ± on iP\ ± ol), (o\ ± o\)\ . 

The operators of interest can be written in this new basis as 



Q 

qT 

Q 



3/2 



20+, 



;b.28) 

;b.29) 
;b.3o) 
;b.si) 
;b.32) 



where the factors of two come from the additional Wick contractions for Q and Qs compared to the 
four flavor operators. Fierz symmetry forbids mixing between the Of and 0~ , and so the matching 
coefficient falls into two 5x5 blocks. We give the results for the Fierz self-conjugate operators: 



£>+(NDR) - £>+(lat) = AcJ (9+(lat) 



(B.33) 
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NDR result 


NDR value 


DRED - NDR 


z+ 


|(5cy + 5c A -c s -cp)-2 


-50.841 -A£ + 2 tad 


3 


z x 


|(9cy + 9c A -C S -Cp) + l 


-47.698 -2£ + 2 tad 


1 

0; 




I{cs + cp) + 1 


-46.398 + 16£ + 2 tad 


5 
3 


Z78 


\{c s + c P -c v -c A ) - 3 


-2.567 + 6£ 


7 
2 


Zs7 


-3 


-3 


3 




i(10c T + llc 5 + llc P -2) 


-46.894 + f£ + 2tad 


7 

6 


z 5 


±(10c T - c s - c P - 2) 


-45.384 - f£ + 2tad 


10 

3 


Z45 


^(2c T - c 5 - c P - 26) 


-2.033 - \i 


2 


Z54 


§(2c T - c s - c P - 2) 


-0.131 - f £ 


4 
3 


z* 


cs -c P = 2(c A - c v ) 


9.644 






Table 11: Matching coefficients, as defined in Eqs. p.33| and p.34| , needed for Fierz self-conjugate 
four-fermion operators. The expressions are given for the NDR scheme and are valid for any local 
lattice operator. Expressions for DRED are obtained by keeping only the c, terms. The numerical 
values are for the NDR scheme with unimproved Wilson lattice fermions, with £ = log(yua), and 
tad = 17.14 being the effect of tadpole improvement. The last column gives the difference between 
the numerical values for DRED and NDR schemes. 



with 



/ 


z + 











\ 




/ 





-22 


2 


16 


16 \ 









Z78 








z* 




-16 


-6 


2 


-24 


8 







Z%7 


z 2 








+ 24 




-4 


2 


-6 


-64 
















z 4 


Z45 




2 


-1 


-13 








V 











Z54 


z 5 ) 




V 


2 


2 


2 








;b.34) 



Our results in the 
t The 



The expressions and numerical results for the Z's are collected in Table [XT 
DRED scheme for the matching coefficients of Oi + , O2+ and agree with those of Ref, 
results of Ref. Q for Z + , Z\ and Z 2 using NDR are, however, wrong, and are corrected by those 
given here. 

The matrix elements needed in the study of B-lifetimes are those of operators Of and O5 . Since 
these are new, we write them in a form which is more directly practical. In the 4-flavor notation 
these operators are 

o% = (s 2 + v 2 )-Us 1 + v 1 -r 1 ), 



{s 1 + v 1 )-Us 2 + v 2 -r 2 ) 



(B.35) 



The 1-loop continuum operators, in the NDR scheme, are 

C^(cont) = [l + A(-46.894 + 2tad + 28/31n(/ia))]0^(lat) 
+A(-2.033 - 4/3 ln(/ia))C+(lat) 

2QS 2 - 2V 1 - 26V 2 - 11V 1 + V 2 + 15.4 1 + 3^ 2 1 



24 

C>+(cont) = [1 + A(-45.384 + 2tad-32/31n( / ua))]C+(lat) 
+A(-0.131 - 16/31n(/ia))0^(lat) 

+^-i(v l - S 1 ) + (V 2 - S 2 ) + V 1 + V 2 ] . 



(B.36) 
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The B parameters are obtained by dividing the above by their VSA forms 

0+(VSA) = (+5/6)[l + AZ P ] 2 (P®P), 

C> 5 + (VSA) = (-1/6)[1 + \Z P f{P®P). (B.37) 

These two operators mix under the renormalization group flow. The 1-loop anomalous dimension 
matrix governing this flow is 

-28/3 4/3 \ 
16/3 32/3 ) 

and the 2-loop contribution is yet to be calculated. Consequently, we state our results for B- 
parameters after the 1-loop matching at scale \i = q*, i.e. without any running in the continuum 
which requires the 2-loop results. 
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